home *** CD-ROM | disk | FTP | other *** search
/ MacFormat 1994 August / August CD.bin / Shareware / Education / numericalmethods Folder / chap_7 / romber.m < prev    next >
Encoding:
Text File  |  1994-06-05  |  876 b   |  35 lines  |  [MATF/MATL]
  1. function [R,quad,err,h] = romber(f,a,b,n,toler)
  2. % [R,quad,err,h] = romber(f,a,b,n,toler)
  3. % Quadrature using the Romberg integration.
  4. % f is the name of the function, input.
  5. % a is the left  endpoint, input.
  6. % b is the right endpoint, input.
  7. % n is the maximum number of rows in the table, input.
  8. % toler is the tolerance, input.
  9. % R is the Romberg table, output.
  10. % quad is the quadrature value, output.
  11. % err is the error estimate, output.
  12. % h is the smallest step size used, output.
  13. m  = 1;
  14. h  = b - a;
  15. err = 1;
  16. j = 0;
  17. R = zeros(4,4);
  18. R(1,1) = h*(feval(f,a) + feval(f,b))/2;
  19. while ((err>toler)&(j<n))|(j<4)
  20.   j = j+1;
  21.   h = h/2;
  22.   s = 0;
  23.   for p = 1:m;
  24.     x = a + h*(2*p-1);
  25.     s = s + feval(f,x);
  26.   end
  27.   R(j+1,1) = R(j,1)/2 + h*s;
  28.   m = 2*m;
  29.   for k=1:j,
  30.     R(j+1,k+1) = R(j+1,k) + (R(j+1,k)-R(j,k))/(4^k-1);
  31.   end
  32.   err = abs(R(j,j)-R(j+1,k+1));
  33. end
  34. quad = R(j+1,j+1);
  35.